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ABSTRACT 

We present collapse simulations of strongly magnetised, turbulent molecular cloud 
cores with masses ranging from 2.6 to 1000 Mq in order to study the influence of 
the initial conditions on the tur bulence-induced disc formation mechanism proposed 



recently by Seifried et al. (2012). We find that Keplerian discs are formed in all cases 
independently of the core mass, the strength of turbulence, or the presence of global 
rotation. The discs appear within a few kyr after the formation of the protostar, are 
50 - 150 AU in size, and have masses between 0.05 and a few 0.1 Mq. During the 
formation of the discs the mass-to-flux ratio stays well below the critical value of 10 
for Keplerian disc formation. Hence, flux-loss alone cannot explain the formation of 
Keplerian discs. The formation of rotationally supported discs at such early phases is 
rather due to the disordered magnetic field structure and due to turbulent motions 
in the surroundings of the discs, two effects lowering the classical magnetic braking 
efficiency. Binary systems occurring in the discs are mainly formed via the disc cap- 
turing mechanism rather than via disc fragmentation, which is largely suppressed by 
the presence of magnetic fields. 
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1 INTRODUCTION 

The formation of protostellar discs under the influence of 



magnetic fields has received a great 


deal of attention for 


about one decade (e.g. Allen et al. 


2003; 


Matsumoto & 


Tomisaka 2004; Machida et al.||2005| |Banerjee k, Pudritz 


2006, |20071 |Price & Bate|2007||Hennebelle fe Fromang|2008 ; 


Hennebelle & Ciardi||2009| |Duffin & Pudritz |2009 


|Com- 


mergon et al.|2010j |Peters et al.|2011 


[Seifried et al. 


2011). 



In these works a highly idealised numerical setup usually 
consisting of an overall core rotation and a magnetic field 
parallel to the rotation axis was used. One common find- 
ing of these simulations is that for an initial setup with a 
magnetic field strength comparable to observations, i.e. a 
mass-to-flux ratio smaller than 5-10 ( e.g. |Falgarone et al. 
2008| |Girart et al.|[2Q09j |Beuther et"aH|2010) , nQ rotation- 
ally supported discs were found. This result is also called 
the "magnetic braking catastrophe" owing to the fact that 
magnetic braking is responsible for the removal of the angu- 
lar momentum which would be necessary to form the discs. 
However, these results seem to be in contrast to recent high- 



resolution observations revealing the existence of Keplerian 
discs around Class protostellar objects 



(e.g. 



2012). 



Tobin et al. 



Trying to circumvent this "magnetic braking catastro- 
phe" by the inclusion of non-ideal MHD effects seems to fail. 



For instance the inclusion of ambipolar diffusion (e.g. Mel- 
lon fe Li|2009| puffin fe Pudritz|2009t does not result in the 
formation of Keplerian discs in the earliest phase of proto- 
stellar evolution. Furthermore, Ohmic dissipation does also 
not allow for large-scale (~ 100 AU) but only for very small 
(~ 10 solar radii) rotationally supported structures (e.g. 



Dapp et al. 2011), unless one assumes an unusually high 



resistivity (iKrasnopolsky et al. 2010). Even the combined 



effects of ambipolar diffusion and Ohmic dissipation cannot 



reduce the magnetic braking efficiency significantly ( |Li et al 
2011). Recently, a third non-ideal MHD effect namely the 
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Hall effect was shown to enable the formation of Keplerian 
discs ( |Krasnopolsky et al.|2011 ). However, this seems to re- 
quire a Hall coefficient about one order of magnitude larger 
than expected under realistic conditions. Furthermore, the 
formation of a disc in these simulations is not due to a re- 
duced magnetic braking efficiency but due to the spin-up of 
the gas in the midplane by the Hall effect as the authors 
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clearly demonstrated by the occurrence of a disc spinning 
in the opposite direction than the surrounding protostellar 



As mentioned before, all the simulations use a highly 
idealised setup of a uniformly rotating core and a magnetic 
field parallel to the rotation axis. Only recently have the 
effects of deviating from such a highly idealised setup been 



included. Hennebelle & Ciardi (2009), Ciardi &; Hennebelle 
(2010) and|Joos et al.| (|2012|) investigated how an inclina- 



tion between the overall magnetic field and the core rotation 
affects the formation of Keplerian discs. The authors found 
that even for a moderate inclination the formation of rota- 
tionally supported discs is possible. However, as Krumholz 
|et al.| ( |2013| pointed out recently, for a realistic distribution 
of magnetic field strengths and misalignment angles, this 
would lead to a fraction of Keplerian discs between 10 - 50% 
only. In an alternative approach, the inclusion of turbulent 
motions in the initial conditions and their effect on the for- 
mation of rotationally supported discs was studied ( Seifried 
et al.||2012| |Santos-Lima et al.||2012a|b| ), again showing the 
possibility of forming Keplerian discs for strong magnetic 
fields. Santos-Lima et ah] ( 2012a|b| ) attribute the formation 
to turbulent reconnect ion occurring in the region of the disc, 
which in turn reduces the magnetic flux and thus the mag- 



netic braking efficiency. In contrast, in Seifried et al. (2012) 
we suggest that the build-up of Keplerian discs is due to the 
turbulent motions in the surroundings of the discs rather 
than due to magnetic flux loss. These motions distort the 
magnetic field and hamper the build-up of a toroidal mag- 
netic field component responsible for angular momentum ex- 
traction, but simultaneously provide a sufficient amount of 
angular momentum necessary for disc formation. 

In this work we investigate this - as we call it - 
turbulence-induced formation mechanism in more detail and 



for a wider range of initial conditions. In |Seifried et al. ( 2012 ) 
we tested this mechanism only for a massive molecular cloud 
of 100 Mq and supersonic turbulence. Since we consider 
the turbulent motions to be the main driver of the forma- 
tion of Keplerian discs, it is particularly interesting to test 
this mechanism for the weaker case of subsonic turbulence. 
As subsonic turbulence is expected to dominate in low-mass 
molecular cloud cores, we extend our analysis also to the 
low-mass range. Moreover, we study if or to what extent the 
turbulence-induced formation mechanism depends on the 
presence of an overall uniform rotation. We show that in 
the presence of turbulence in essentially any case Keplerian 
discs are formed by the aforementioned mechanism. 

The paper is organised as follows. In Section [2] we 
present the details of the initial conditions of the simulations 
and briefly describe the numerical techniques. The result 
of the simulations are presented are in Section [3] First, we 
present the details of the Keplerian discs formed and anal- 
yse possible reasons for their formation. Next, we consider 
the time evolution of global disc properties before studying 
the fragmentation properties of the discs in more detail. In 
Section [4] the results are discussed in a broader context and 
are compared to related numerical and observational work 
before we summarise our main findings in Section [5] 



2 INITIAL CONDITIONS 

We now describe the setup of the different simulations. We 
performed several simulations with core masses ranging from 
to 2.6 Mq to 1000 M^ (see Table [Tl for an overview over 



A 



the different models and their parameters) . All cores are em- 
bedded in a cubic box with a length of approximately three 
times the cores diameter to avoid corruption of the results 
by boundary effects. The ambient medium has a density 
a factor of 100 lower than the density at the edge of the 
cores, hence its dynamical influence can be considered neg- 
ligible. The runs 2.6-NoRot-M2 and 2.6-Rot-M2 have a core 
mass of 2.6 Mq and a gas temperature of 15 Krj Here we 
note that the first part in the name of each run denotes the 
mass-to- flux ratio, the second whether a uniform rotation is 
present or not, and the third the (approximate) core mass 
in Mq. The remaining runs with higher core masses have a 
initial temperature of 20 K, thus somewhat higher than in 
the low-mass cores (see e.g. |Ragan et al.||2012) |Launhardt 
|et al.||2013 ). 

In all runs we set up a Kolmogorov-type turbulence field 
so that 3D-turbulent rms Mach number (M rms ) ranges from 
subsonic turbulence (M rms ~ 0.74) in the low-mass runs 
to strongly supersonic turbulence in the high-mass cases 
(M rms ~ 2.5 — 5.4). These values represent the typical level 
of turbulence observed in molecular cloud cores (e.g. Caselli 



|fe Myers|[l995| |Andre et al.||2007| but see also the review 
of Ward-Thompson et al.||2007 and references therein) and 
result in a ratio of the turbulent energy to gravitational en- 
ergy /3turb of about 0.08 to 0.12 (see Table[T]). Here the gravi- 
tational energy is determined from the analytical expression 
for the given density profile (see further down J] For the runs 
with 100 Mq we performed several simulations in order to 
test to what extent the results depend on the random reali- 
sations of the turbulence field (seed A,B, and C). We found 
that the results do not change qualitatively if the turbulence 
field is varied (but see |Seifried et al.|2012 for a detailed dis- 
cussion). In some of the runs in addition we superimpose 
a rigid rotation of the core around the z-axis on the tur- 
bulent velocity field. For the sake of simplicity the rotation 
frequency in these runs is chosen such that the rotation en- 
ergy equals the turbulent kinetic energy, i.e. f3 vo t — Aurb- All 
cores are threaded by a magnetic field parallel to the z-axis. 
The strength of the field declines with increasing cylindrical 
radius R as 



B z oc V^OR) 



midplane 5 



(1) 



where p(i?) m id P iane is the density in the midplane at the 
radius R. This guarantees that /3 p i asma = P/(B 2 /Sir) is con- 
stant in the equatorial plane, where P is the thermal pres- 
sure. Furthermore, in order to guarantee VB = 0, B z does 
not vary along the z-axis. The strength of the magnetic field 
is chosen such that the (normalised) mass-to-flux ratio fi of 



1 The ambient gas has a 100 times higher gas temperature to 
assure pressure equilibrium at the edge of the core. 

2 We note that in our previous work a value for (3 of 0.04 was 
given. This value resulted from the calculation of the gravitational 
energy from the simulation data not taking into account a sys- 
tematic offset in the numerical gravitational potential compared 
to the analytical value, which results in the somewhat lower value 
for p. 
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Table 1. Initial conditions of the performed simulations showing the mass, radius, normalised mass-to-flux ratio fi, whether uniform 
rotation is present or not, the angular frequency, the turbulent energy content normalised to the gravitational energy, the seed of the 
random realisation of the turbulence field, the power spectrum index p of the turbulence spectrum, the rms Mach number, and simulated 
physical times of the protostellar discs. 



Run 


T^core 


r core 


/i 


rotation 


n 


Pturb 


turbulence 


P 


M rms 






[M ] 


[pc] 






[10- 13 s- 1 ] 




seed 






[kyr] 


2.6-NoRot-M2 


2.6 


0.0485 


2.6 


No 





0.087 


A 


5/3 


0.74 


15 


2.6-Rot-M2 


2.6 


0.0485 


2.6 


Yes 


2.20 


0.087 


A 


5/3 


0.74 


15 


2.6-NoRot-M100 


100 


0.125 


2.6 


No 





0.084 


A 


5/3 


2.5 


15 


2.6-Rot-M100 


100 


0.125 


2.6 


Yes 


3.16 


0.084 


A 


5/3 


2.5 


15 


2.6-Rot-M100-B 


100 


0.125 


2.6 


Yes 


3.16 


0.084 


B 


5/3 


2.5 


15 


2.6-Rot-M100-C 


100 


0.125 


2.6 


Yes 


3.16 


0.084 


C 


5/3 


2.5 


15 


2.6-Rot-M100-p2 


100 


0.125 


2.6 


Yes 


3.16 


0.084 


A 


2 


2.5 


15 


2.6-NoRot-M300 


300 


0.125 


2.6 


No 





0.12 


A 


5/3 


5.0 


10 


2.6-Rot-M1000 


1000 


0.375 


2.6 


Yes 


1.90 


0.081 


A 


5/3 


5.4 


10 



the cores is 2.6. We emphasise that throughout the paper 
the mass-to-flux ratio is given in units of the critical mass - 
to-flux ratio /i C rit = 0.13 /y/G ( |Mouschovias k, Spitzer|l976 ), 
where G is the gravitational constant. 

The density profile varies between the different runs. In 
the runs 2.6-NoRot-M2 and 2.6-Rot-M2 the density profile 
follows that of a Bonnor-Ebert (BE) sphere up to a radius 
of 6.9ro (6.49ro is the critical radius for a BE sphere). The 
physical radius of the core is 10 000 AU, and the density of 
the BE-profile is scaled up in order to contain 2.6 Mq in 
the core. Therefore, the core is slightly gravitational unsta- 
ble, which is reasonable since additional effects like rotation, 
turbulence and magnetic fields can counteract gravity. The 
runs with 100 Mq and 300 Mq have a diameter of 0.125 pc 
and are highly gravitationally unstable (56 and 290 Jeans 
masses, respectively). The density of the cores declines out- 
wards ajf] 



p(r) oc r 



(2) 



Since the cores are highly gravitationally unstable, this ex- 
plains the higher turbulent Mach numbers compared to the 
runs with 2.6 Mq despite comparable values of /3 tU rb- Run 
2.6-Rot-M1000 has a mass of 1000 Mq corresponding to 
about 340 Jeans masses contained in a core with a radius 
of 0.375 pc having qualitatively the same density profile as 
described in Eq. [2] 

The cooling routine applied in all runs takes into ac- 
count dust cooling, molecular line cooling and the effects of 
optically thick gas (Banerjee & Pudritz 2006). We introduce 
sink particl es above a density thr eshold of p C rit = 1.14- 10 -10 

for details) . The maximum 



Federrath et al. 



2010 



g cm (see 

spatial resolution is set to 1.2 AU in all simulations. The re- 
finement criterion used guarantees that the Jeans length is 
resolved everywhere with at least 8 grid cells. 



3 To avoid unphysically high densities in the interior of the core, 
we cut off the r -1-5 profile at a radius of 0.0125 pc. Within this 
radius the density distribution follows a parabola p(r) oc [1 — 
(r/ro) 2 ]- 



3 RESULTS 

In what follows we focus on the effects of varying core masses 
and initial rotation on the formation of early- type discs. For 
this reason we present details of the runs 2.6-NoRot-M2, 2.6- 
Rot-M2, 2.6-NoRot-M100, 2.6-Rot-M100, 2.6-NoRot-M300 
and 2.6-Rot-M1000 only. The remaining runs with core 
masses of 100 Mq listed i n Table [l] were discussed in de- 
tail in |Seifried et al. ( 2012| . We again emphasise that we 
did not find any qualitative differences for different realisa- 
tions of the turbulence field. We therefore believe that also 
in the runs with lower and higher masses a change in the 
turbulence field would not significantly alter the results (see 
also Joos et aT]|2013| . 



3.1 Disc formation 

3.1.1 Column density structure 

Our first step is to compute the column density structure of 
the discs formed in the simulations and to analyse how they 
evolve over time. We note that the column density is com- 
puted for the discs for the face on view. For the calculation 
of the disc's angular momentum we only take into account 
gas with a density larger than 5 • 10 -13 g cm -3 . Next, the 
angular momentum is calculated with respect to the centre- 
of-mass of the disc. The threshold density of 5- 10 -13 g cm -3 
is physically motivated by the fact that gas gets optically 
thick in this range (see also the discussion in Section 4.3). 



We found that when using this threshold the radius of the 
disc also reasonably well agrees with the rather sharp drop- 
off in the column density marking the outer edge of the disc 
(see figures below). Furthermore, making use of a density 
threshold has the advantage that gas from the outflow cav- 
ity above and below the disc, which usually has significantly 
lower densities, is excluded. 

We plot the column density structure of the six runs 
mentioned in the beginning of this section for two differ- 
ent times, namely half way from the formation of the first 
sink particle towards the end of the simulation (7.5 and 5 
kyr, respectively) and directly at the end of each simulation 
(15 and 10 kyr, respectively). In addition, we also show the 
line-of-sight averaged velocity in the discs. First we show 
the results for the discs of the runs with 2.6 Mq cores in 
Fig. [I] The viewing direction was adjusted to the direction 
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Figure 1. Column density in logarithmic scaling for the discs 
in run 2.6-NoRot-M2 (top) and 2.6-Rot-M2 (bottom) for 7.5 kyr 
(left) and 15 kyr (right) after the formation of the first sink par- 
ticle. The viewing direction was adjusted in order to see the discs 
face on. The green arrows show the mass-weighted mean velocity 
and the black dots the projected positions of the sink particles. 
The figures are 200 AU in size. 



of the angular momentum vector such that the discs are 
seen face on. Due to the rather weak turbulent motions and 
the marginally gravitationally unstable configuration of the 
cores, in both runs no fragmentation has occurred and only 
one sink particle with an associated disc has formed so far. 
Furthermore, it can be seen that in both runs a rotationally 
supported structure has built up after 7.5 kyr and that the 
disc in run 2.6-Rot-M2 is somewhat larger that that in run 
2.6-NoRot-M2. 

In Fig. [2] we plot the structure of the discs in run 2.6- 
NoRot-MlOO and run 2.6-Rot-M100. Due to the supersonic 
turbulent motions and the highly gravitationally unstable 
configuration, the cores in both runs have fragmented and 
formed 5 sinks by the end of each run. However, only one disc 



has formed so far in both runs (see Section 3.5 for the de- 
tails of the fragmentation process). In run 2.6-Rot-M100 the 
three most massive sinks are all contained within one mas- 
sive disc. The two remaining sink particles are a few 1000 
AU off from the centre of this disc and do both not contain 
any associated disc, i.e. no gas with densities above 5 • 10 -13 
g cm -3 has been found in their surroundings. We attribute 
this to the fact that both sink particles have been ejected 
from the massive disc before and are thus moving with a 
high relative speed compared to the surrounding medium, 
which prevents the build-up of discs. Similar holds for run 
2.6-NoRot-M100 where one sink has been ejected and one 
sink has been formed shortly before the end of the simulation 
so that nothing can be said about an associated disc, yet. 
The remaining three sink particles are again grouped in one 
disc. Like in the runs with 2.6 Mq, in both runs the discs 
show significant rotationally support already at 7.5 kyr. Fur- 
thermore, it can be inferred from Fig. [2] that in both runs 




Figure 2. Same as in Fig. [T] but for run 2.6-NoRot-M100 (top) 
and run 2.6-Rot-M100 (bottom). 



two of the sink particles create a very tight binary with a 
separation of about 1 AU (not resolved in the figure). The 
third sink in the disc is at a larger distance of about 10 AU. 

In Fig. [3] we consider run 2.6-NoRot-M300, which has 
a highly gravitationally unstable configuration and strongly 
supersonic turbulence with a rms-Mach number of about 5. 
Hence, it is not surprising that run 2.6-NoRot-M300 shows 
heavy fragmentation with more than 50 sink particles hav- 
ing been formed so far. Four discs can be found, of which we 
plot the two most massive ones at a time of 5 kyr and 10 kyr. 
As can be seen, the first disc (top panel of Fig. |3| has been 
subject to massive fragmentation. It contains 8 sink parti- 
cles and shows several very prominent spiral arms. 5 of the 
8 sinks form a very tight and strongly bound system with a 
size of a few AU. Despite the strongly fragmented structure 
the disc still retains a well-defined rotation structure. By the 
end of the simulation the second disc is by far smaller than 
the first one although also in this case a rotationally sup- 
ported structure is recognisable up to a radius of about 25 
AU. Furthermore, it shows significantly less fragmentation 
than the first disc and only contains 2 sink particles. 

Finally, in Fig. [4] we show the two most massive discs 
formed in run 2.6-Rot-M1000. As in run 2.6-NoRot-M300 
the 1000 Mq core has fragmented heavily forming 23 sink 
particles by the end of the simulation. The somewhat lower 
degree of fragmentation (about a factor of 2) is most likely 
a consequence of the much more compact configuration 
with a rsj 4 times higher central density in run 2.6-NoRot- 
M300. The Jeans masses given in Section |2] (340 m run z.o- 
Rot-MlOOO compared to 290 in run 2.6-NoRot-M300) are 
only globally averaged quantities, i.e order-of-magnitude es- 
timates, which can therefore not be taken as an absolute 
indication of the expected degree of fragmentation. As can 
be seen, the most massive disc (top panel of Fig. |4j) is signifi- 
cantly less fragmented than the most massive disc in run 2.6- 
NoRot-M300 (top panel of Fig. [3}, containing only 5 sinks. 
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Table 2. Simulation properties at the end of each run: total num- 
ber of sinks and discs, diameter, mass, and number of fragments 



Figure 3. Same as in Fig. ^ but for the two most massive discs 
(top and bottom row) in run 2.6-NoRot-M300 and times of 5 kyr 
and 10 kyr after the formation of the first sink particle. 




Figure 4. Same as in Fig. [3] but for the two most massive discs 
(top and bottom row) in run 2.6-Rot-M1000. 



A ring with high column density appears at a distance of 
20 - 30 AU from the disc centre which might probably frag- 
ment in the future. The second disc (bottom panel of Fig.[4| 
has not yet formed around 5 kyr and only at the end of the 
simulation a rotating structure up to a radius of about 25 
AU becomes visible. 

Two important results can be inferred from the Figs.[l]- 
[4] First, for all discs the typical radius, i.e. the radius where 
the column density reveals a sharp drop-off, is in general not 
larger than ^50-75 AU (except maybe the first disc in run 



in the most massive 


disc and the mass- 


to-flux ratio in 


its vicinity. 


Run 


N s inks 


Ndiscs 





^disc 


Nfrag 










[AU] 


[M ] 






2.6-NoRot-M2 


1 


1 


67 


0.0585 


1 


6.2 


2.6-Rot-M2 


1 


1 


89 


0.0607 


1 


6.6 


2.6-NoRot-M100 


5 


1 


124 


0.175 


3 


12.9 


2.6-Rot-M100 


5 


1 


124 


0.452 


3 


14.4 


2.6-Rot-M100-B 


36 


4 


71 


0.062 


4 


11.5 


2.6-Rot-M100-C 


5 


1 


111 


0.263 


3 


9.3 


2.6-Rot-M100-p2 


5 


1 


193 


0.63 


4 


14.6 


2.6-NoRot-M300 


57 


4 


132 


0.048 


8 


9.2 


2.6-Rot-M1000 


23 


4 


109 


0.40 


5 


9.7 



2.6-NoRot-M300). This is significantly smaller than the size 
of the rotationally supported discs we find in similar runs 



without turbulence (Seifried et al. 2011). In this work we 



find Keplerian discs with radii up to a few 100 AU although 
for an initial mass-flux-ratio larger by about a factor of 5 or 
more than that in the simulations presented here. A second 
important result inferred from the column density plots is 
the fact that all discs seem to build up a rotationally sup- 
ported structure already during the earliest phase of their 
evolution. Indeed, we checked that in all runs rotationally 
supported structures become visible already around 5 kyr 
after the formation of the first sink particle, thus at an even 
earlier time than that shown in Figs. [I] and [2] 

Some global properties of the discs are listed in Table [2] 
For the mass of the disc the density threshold criterion men- 
tioned in the beginning is used. The radius of the disc is 
defined as the distance from the centre-of-mass where the 
azimuthally averaged column density drops below a value 
of 1 • 10 25 cm -2 . Although this choice is somewhat arbi- 
trary, we found it to be reasonable from visually inspecting 
the column density plots. Furthermore, we only denote the 
accumulation of dense gas around a sink particle as a disc 
in case a clear rotationally supported structure larger than 
about 30 AU is recognisable. 



3.1.2 Velocity structure 

Next, we consider the velocity structure in the discs in more 
detail. For this purpose we calculate the rotation velocity 
and the radial velocity v ra d for each cell in the disc. We note 
that the calculation is done in the rest-frame of the discs, 
i.e. the z-axis is defined by the disc angular momentum vec- 
tor and the velocity is corrected by the velocity of the disc's 
centre-of-mass before calculating and v Ta d • In contrast to 
the column density plots presented before we restrict the 
consideration to the final snapshot at the end of each sim- 
ulation. However, we emphasise that when considering the 
velocity structure at earlier times, the results do not change 
qualitatively as could already be inferred from the velocity 
vectors shown in the column density plots. Furthermore, we 
note that for radii larger than ~ 50 AU (depending on the 
actual simulation), where no cells with gas densities above 
5 • 10 -13 g cm -3 are found, we adopt a simple geometrical 
criterion considering all the gas within a height of 20 AU 
above/below the midplane defined by the disc. Hence, we 
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Figure 5. Radial dependence of the rotation (red) and radial 
velocity (green) for the four discs of the runs 2.6-NoRot-M2, 
2.6-Rot-M2, 2.6-NoRot-M100 and 2.6-Rot-M100 (from top left 
to bottom right) at the end of each simulation. The black solid 
line shows the Keplerian velocity ^kep- The regions below 4 AU 
are affected by resolution effects, therefore they are shaded grey 
to guide the reader's eye. Note the different scaling of the y-axis. 
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Figure 6. Same as in Fig. [5] but for the two most massive discs 
of the runs 2.6-NoRot-M300 (top) and 2.6-Rot-M1000 (bottom). 



are able to analyse the velocity structure also outside the 
disc denned by the density threshold. In order to get an im- 
pression of the scatter of and v ra d, the values are not 
averaged azimuthally. The Keplerian velocity (GM/r) 1 ^ 2 is 
calculated by taking into account the mass of gas and all 
sink particles within a sphere of radius r around the centre 
of the disc. 

In Fig. [5] we plot the radial dependence of and v Ta d 
for the four discs in the runs 2.6-NoRot-M2, 2.6-Rot-M2, 
2.6-NoRot-M100 and 2.6-Rot-M100 and in Fig.[6]for the two 
most massive discs in the runs 2.6-NoRot-M300 and 2.6-Rot- 
M1000. It can be seen that within the typical disc radius 
inferred from the column density plots (< 50 AU) the discs 
in general rotate with the Keplerian velocity (black lines). 
The radial velocity on the other hand scatters around zero 
as expected for a Keplerian disc and is in general signifi- 
cantly smaller than the rotation velocity. Interestingly, for 
some of the runs the Keplerian velocity structure extends 
up to ~ 100 AU and thus even beyond the disc radius in- 



ferred from the column density plots. This result is in strong 
contrast to a series of star formation simulations in which 
turbulence is neglected ( e.g.|Allen et al.|2003 |Price &; Bate 



2007 1 IMellon fe Li|[2008] |Hennebelle Fromang||2008| Puf- 



fin fe Pudritz|2009| jSeifried et al.|2011| ). In general, in these 
simulations it was found that for mass-to-flux ratios below 5 
- 10, i.e. for moderately strong magnetic fields, no rotation- 
ally supported discs are formed. Again, we stress that these 
results are due to the particular choice of the initial condi- 
tions namely the omission of turbulence. On the other hand, 
if turbulence is taken into account, recently the build-up of 
early-type Keplerian discs was reported even for strong mag- 
netic fields ( |Seifried et al.|2012|[Santos-Lima et al.|2012a|b| . 

As demonstrated in this section, the build-up of Kep- 
lerian discs does not depend on the mass of the prestellar 
core - and hence also not on the strength of the turbulence, 
which ranges from subsonic motions for the low-mass cores 
to highly supersonic motions for the most massive cores. Fur- 
thermore, for the 2.6 Mq and 100 Mq cores it can be seen 
that the lack of an overall core rotation does not significantly 
affect the build-up of Keplerian discs since turbulence always 
produces local shear flows which carry a sufficient amount 
of angular momentum to build up the discs. For both cases, 
the rotating and the non-rotating, Keplerian discs build up 
in the low-mass and high- mass case. However, a close in- 
spection Fig. [l] and [2] shows that in the case of an overall 
core rotation the discs seem to be slightly larger than for 
the non-rotation case. This statement, however, has to be 
taken with caution since we obviously require a larger num- 
ber of discs in order to obtain statistically significant results. 
Nevertheless, our result of slightly larger discs for uniformly 
rotating cores agrees with what one would naively expect. 
To summarise, the mechanism of turbulence-induced 



disc formation as proposed recently in Seifried et al. (2012) 



works for a wide range of core masses, turbulence strengths, 
and even for cores without uniform initial rotation. We again 
emphasise that the Keplerian disc structure seen in Figs. [5] 
and [6] at 15 kyr is present at much earlier times already, i.e. 
as early as 5 kyr. Assuming rotation velocities of a few km 
s _1 at a radius of 100 AU, this time corresponds to about one 
to a few rotation periods of the disc - and correspondingly 
more at smaller radii- indicating that Keplerian discs build 
up in the very beginning. In the following we will explore this 
mechanism and the differences to the non-turbulent case in 
more detail. 



3.2 Flux loss 

The extensive work on disc formation under the influence of 
magnetic fields has revealed a critical value for the mass-to- 
flux ratio of 5 - 10 below which the formation of Keplerian 
discs is suppressed in collapsing cloud cores without tur- 
bulent motions (e.g. |Allen et al 



2003 



Price & Bate 



2007 



Mell on fe L1U2OO8I [Hennebelle &: Fr omang 2008; Puffin fe| 



Pudritz|2009[ |Seifried et al.||2011| ). In the following we test 
to what extent the formation of Keplerian discs in our sim- 
ulations can be explained by the loss of magnetic flux in the 
discs' vicinity. For this purpose we calculate the mass-to-flux 
ratio fi in spheres of radius r = 500 AU around the centre- 
of-mass of the discs formed in the four runs 2.6-NoRot-M2, 
2.6-Rot-M2, 2.6-NoRot-M100 and 2.6-Rot-M100. We note 
that we have chosen a radius of 500 AU for a particular rea- 
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Figure 7. Top: Mass-to-flux ratio for spheres with a radius of 
500 AU around the discs of the runs 2.6-NoRot-M2, 2.6-Rot-M2, 
2.6-NoRot-M100 and 2.6-Rot-M100. The straight, black, dashed 
line indicates the threshold for /i below which Keplerian disc for- 
mation is suppressed. The arrow indicates the time when a disc 
with significant rotational support is observed for the first time. 
Bottom: Inclination of the disc angular momentum to the mean 
magnetic field for the same spheres as in the top panel. 



son: By calculating the gas torque and the magnetic torque 
m |Seifried et aL] ( |2Qll| ) we have shown that a large fraction 
of the angular momentum is removed already before the gas 
falls onto the disc. The analysis showed that up to a distance 
of 500 - 1000 AU from the centre of the disc the (negative) 
magnetic torque balances the (positive) inwards angular mo- 
mentum transport by the gas (see bottom panel of Fig. 11 
in |Seifried et al.||2011| ). Hence, this result clearly illustrates 
the importance to check what happens on scales larger than 
100 AU. 

In the top panel of Fig. [7] we plot \i in the r = 500 
AU-sphere. It can be seen that fi increases over time reach- 
ing values > 10 for the 100 Mq runs. This indicates that 
there is a loss of magnetic flux in the vicinity of the discs 
over time by a factor of a few. However, we emphasise that 
a part of the increase of \i is most likely due to accretion 
along the magnetic field lines - a fact which cannot be cov- 
ered by a simple geometrical calculation of the mass-to-flux 
ratio and which is often ignored in literature. Hence, it is 
likely that the influence of flux loss e.g. by magnetic recon- 
nection (Santos-Lima et al. 2012a b) on the increase of \i 
and thus the formation of Keplerian discs is somewhat over- 
estimated. 

Despite the steady increase of \i the observed flux loss is 
still too small to explain the build-up of Keplerian discs on 
its own. For the 100 Mq runs \i exceeds 10 only after 10 kyr, 
i.e. well after the Keplerian disc has built up, and for the 2.6 



Mq runs it stays below 10 for the entire 15 kyr. Since non- 
turbulent simulations, however, suggest that a value of \i > 
10 is necessary for Keplerian discs to form, clearly the loss 




Figure 8. Edge-on view of the disc in run 2.6-Rot-M100. The 
column density is overplotted with the magnetic field lines (white) 
and the sink particles (green dots). The figure is 1000 AU in size 
thus showing the same region analysed in Fig. [7] 



of magnetic flux alone cannot explain the observed build- 
up of Keplerian discs at t < 5 kyr. We note that for the 
remaining runs not considered here the time evolution of fi 
is qualitatively and quantitatively very similar to the runs 
2.6-NoRot-M100 and 2.6-Rot-M100, i.e. fi reaches values > 
10 at the end of the simulation (see Table 



3.3 Magnetic field structure 



Recently, [Hennebelle Ciardi (2009), Ciardi & Hennebelle 



(2010) and Joos et al. (2012) suggested that an inclination 



of the magnetic field with respect to the rotation axis of the 
core decreases the magnetic braking efficiency. In order to 
test this statement, we inspected the magnetic field struc- 
ture in the vicinity of our discs. For this purpose, in the bot- 
tom panel of Fig. [7] we plot the inclination a of the disc an- 
gular momentum to the average magnetic field in the sphere 
of r = 500 AU. The inclination a is found to be relatively 
large except for run 2.6-Rot-M2. We note that also for the 
remaining runs there is in general a significant inclination 
between the disc angular momentum and the mean magnetic 
field. Interestingly, a seems to be somewhat smaller for the 
runs with overall rotation. This can be easily understood 
given the fact that for these runs the disc angular momen- 
tum at later times is much closer to the z-axis (see left panel 
of Fig. [9} and thus also closer to the overall direction of the 
magnetic field in the core. 

At first sight, the relatively large inclination between 
the angular momentum of the disc and the mean mag- 
netic field support the findings of |Hennebelle &; C iardi 
( [2009t , [Ciardi fe Hennebelie] ( pOlOf and |Joos et al.| fl2012| ) 
showing that inclined magnetic field lines aid the formation 
of Keplerian discs even for high magnetic field strengths. 
To investigate those arguments in more detail, we inspect 
the structure of the magnetic field lines in the vicinity of 
the disc in run 2.6-Rot-M100 in Fig. [§] Although a pre- 
ferred direction of the magnetic field is recognisable, there 
is a significant distortion of the magnetic field due to the 
turbulence on scales of a few 100 AU. For such a distorted 
magnetic field, however, the definition of a mean field direc- 
tion, as used by |Hennebelle &; C iardi (2009|, [Ciardi &; Hen-| 



8 D. Seifried et al 



nebelle (2010) and Joos et al. (2012), and its inclination to 



the disc orientation can be considered questionable. Hence, 
we strongly suggest that the large inclination between the 
disc orientation and the mean magnetic field measured in 
the bottom panel of Fig. [7] is most likely not the only reason 
for the build-up of Keplerian discs. 

The distortion of the magnetic field in the surround- 
ings of the discs itself can explain the reduced magnetic 
braking efficiency and the subsequent build-up of Keplerian 
discs since the disordered field structure prevents an effec- 
tive coupling of the inner parts close to the disc to the outer 
parts. A magnetic field line structure, however, which fans 
out from a typical disc radius Ro to a radius R somewhere 
outside, would significantly increas es the magnetic br aking 
efficiency by a factor of (R/Rq) 4 (Mouschovias 1985). For 



a non-turbulent simulation usually such a configuration is 
present, see e.g. Fig. 13 in |Hennebelle &; Fromang (2008} or 



Fig. 9 in Seifried et al. (2011). In our simulations, however, 



only a strongly disordered field structure is present due to 
the turbulent motions. This is what we consider to be the 
main reason for the reduced magnetic braking efficiency. We 
emphasise, however, that more generally one could regard 
the reduced magnetic braking efficiency due to misaligned 



magnetic fiel ds (jHennebelle &; Ciardi 2009 Ciardi &; Hen 
nebelle 2010} Joos et al. 2012) as being closely related to 



our case when interpreting the misalignment as some kind 
of disorder. 

Moreover, the magnetic braking mechanism as de- 



scribed by Mouschovias & Paleologou (1980) relies on the 



build-up of a toroidal magnetic field due to coherent rota- 
tional motions. In the vicinity of the discs presented here 
such coherent rotation structures are barely present. In- 
stead, mass accretion occurs along distinct channels (shear 
flows) created locally by the turbulent motions (see Figs. [I] 
- [4} . Since we have shown in Seifri ed et al.| ( |2011 ) that a 
significant fraction of the angular momentum is removed by 
magnetic braking already before the gas falls onto the disc, 
the lack of a coherent rotation structure in the surroundings 
of the discs further reduces the magnetic braking efficiency. 
Simultaneously, the shear flows can provide enough angu- 
lar momentum to build up the discs as already pointed out 



Seifried et al. (2012) 



3.4 Time evolution 

For the analysis of the time evolution of the discs we re- 
strict ourselves to the runs 2.6-NoRot-M2, 2.6-Rot-M2, 2.6- 
NoRot-MlOO and 2.6-Rot-M100. There are some interesting 
differences between low-mass and high-mass cores as well as 
rotating and non-rotating cores, which will be demonstrated 
in the following. We remind the reader that we define the 
discs by all the gas with a density larger than 5 • 10 -13 g 
cm -3 around the corresponding sink particle. 

First, we consider the time evolution of the orientation 
of the discs in the left panel of Fig. [9] For this purpose we 
calculate the angular momentum of the discs and show its 
deviation from the z-axis, i.e. the polar angle 0. We have cho- 
sen the z-axis as it corresponds to the direction of the initial 
magnetic field. As can be seen there is a significant difference 
between the runs with an overall core rotation compared to 
that without any rotation as well as between the low- and 
high-mass runs. In the beginning the four discs have an an- 



gular momentum well off from the z-axis. However, for the 
runs 2.6-Rot-M2 and 2.6-Rot-M100 the orientation of the 
disc angular momentum approaches the z-axis (0 ~ 0) after 
some time, for run 2.6-Rot-M2 significantly earlier than for 
run 2.6-Rot-M100. Moreover, one can see that this change 
in the direction of the angular momentum happens rather 
quickly - within about 1000 yr. After the disc orientation has 
approached the z-axis it stay more or less constant, which 
is particularly interesting for run 2.6-Rot-M100 since here 
strongly supersonic turbulent motions are present. Such an 
alignment of the rotation and magnetic field direction was 



first observed by Machida et al. ( 2006 ) 



Interestingly, for run 2.6-NoRot-M2, where no overall 
rotation but only turbulent motions are present, the disc 
orientation remains remarkable constant and also the az- 
imuthal angle <j> does not change significantly over time. We 
note that we decided not to show the azimuthal angle in 
Fig. [5] since for the cases where is close to zero the varia- 
tion in <j) can be very large - although the disc orientation it- 
self does hardly change in this case - which would make the 
figure basically unreadable. In contrast to run 2.6-NoRot- 
M2, run 2.6-NoRot-M100 shows large variations in the disc 
orientation, which we attribute to the very turbulent envi- 
ronment they form in. Again the changes in the disc's orien- 
tation occur rapidly within about 1000 yr or even less. This 
points to very distinct accretion events of blobs of gas onto 
the disc, which results in an almost instantaneous change of 
the angular momentum. 

In run 2.6-Rot-M100 the mass of the disc is a factor of 
a few larger than that of the disc in run 2.6-NoRot-M100. 
This could be attributed to the closer orientation of the disc 
in run 2.6-Rot-M100 to the z-axis thus profiting from the 
large reservoir of angular momentum due to the overall ro- 
tation. For the runs with 2.6 Mq, however, both discs have 
rather similar masses of about 0.05 Mq. Furthermore, these 
discs do not grow over time but keep their mass whereas 
the discs in the high-mass core runs show a clear increase 
of mass over time despite distinct events of significant mass 
loss. Moreover, a close comparison between the disc orien- 
tation and mass for run 2.6-NoRot-M100 shows that the 
rapid change of the disc's orientation after 4 kyr and 6 kyr 
is correlated to a strong increase of the disc mass. This con- 
firms our assumption of distinct accretion events onto the 
disc being responsible for the rapid changes. We note that 
at the end of each simulation the disc masses are well be- 
low the total sink particle masses in the discs, which range 
from 0.25 M in run 2.6-Rot-M2 to 6.3 M in run 2.6- 
NoRot-MlOO (right panel of Fig. [9). We find that phases of 
significant decrease in disc mass (see middle panel of Fig. [9]) 
are usually accompanied with increased sink accretion rates. 
This can be explained by the preceding accretion of gas onto 
the disc with a small or even negative angular momentum 
(with respect to the disc angular momentum). Hence, the 
specific angular momentum of the disc is lowered and so 
is the rotational support, which in turn results in a subse- 
quent increase of the sink accretion rate until a rotationally 
supported equilibrium in the disc is established again. 

We emphasise that the discs in the runs not discussed 
above show qualitatively very similar properties with typical 
masses of a few 0.1 Mq and partly strong fluctuations (see 
Table [2]) and a disc angular momentum well off from the 
z-axis (and thus from the overall magnetic field). However, 
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Figure 9. Left: Time evolution of the disc orientation as measured by the angle between the z-axis and the disc angular momentum. 
For the runs with initial uniform rotation the disc's angular momentum clearly approaches the overall angular momentum over time. 
Middle: Time evolution of the disc masses. Right: Time evolution of the total mass of all sinks contained in the disc. 



as for the runs 2.6-Rot-M2 and 2.6-Rot-M100 in all other 
runs with an overall core rotation the orientation of the disc 
approaches the z-axis over time as seen already before. The 
only exemption is run 2.6-Rot-M100-C, where the disc stays 
misaligned over the entire 15 kyr. 

To summarise, we find that there are some differences 
in the time evolution of the discs between the runs with 
and without core rotation as well as between high- and low- 
mass cores. Nevertheless, as demonstrated in Section pTT) the 
velocity structure of the individual discs is barely affected. 

3.5 Fragmentation behaviour 

Recently, it was suggested that a large fraction of low-mass 
protostar and brown dwarfs form via gravitational fragmen- 
tation of protostellar discs (see e.g. Stamatellos et aL|2011[ 
and references therein). In contrast, Bate (2009) and Offner 
et al. (2009) find that the degree of disc fragmentation is 



reduced when including radiative feedback. In these simu- 
lations low-mass stars mainly form via the collapse of dis- 
tinct gravitational cores or condensations (see also Padoan 



& Nordlund 2004 Hennebelle & Chabrier 2008). At first 



sight the multiple sink particle systems in the discs in our 
high-mass core runs support the classical disc fragmentation 
picture ( |Toomre] fl964). However, a second mechanism to 
form such a system is the so-called disc-assisted capture (see 
e.g. Section 5.4 in Zinnecker & Yorke 2007, and references 
therein) where the cross section for capturing a bypassing 
protostar is increased significantly due to the presence of an 
extended disc. In the following we examine the exact mech- 
anism how the multiple protostellar systems are formed in 
our simulations. 

In Fig |10| we show the relative distance of each sink 
particle to the sink particle formed first for run 2.6-Rot- 
M100, i.e. 



4(t) = |ri(t)-n(t)| : 



(3) 



where r*i(t) is the position of the i-th sink particle at time 
t and r*i (t) that of the first sink particle. As can be seen, 
by the time the sink particles form they are more than 100 
AU away from the first sink particle and the centre of its 
associated disc. Hence, these sink particles do not form by 
fragmentation of the already existing disc but clearly out- 
side it during the collapse of distinct, gravitational unstable 
regions. Only during their further evolution do the secon- 
darily formed sinks approach the first sink particle and are 




t/kyr 



Figure 10. Distance of the different sink particles to the sink 
particle formed first (see Eq. [3} for the simulation 2.6-Rot-M100. 



then captured in the disc. From Fig. [l0| it can also be seen 
that two of the sinks get ejected from the disc due to close 
three-particle interactions and move away relatively quickly. 
These ejected sinks usually do not build up an associated 
disc around them as already pointed out in Section [3T] 

At first sight it seems interesting that the sink particles 
in run 2.6-Rot-M100 do not form via disc fragmentation but 
rather due to the fragmentation of the entire core. For this 



reason we consider the Toomre parameter (Toomre 1964) 
describing the stability of the disc: 

KC S 



Q = 



tt^G 



(4) 



with the epicyclic frequency k, sound speed c s , surface den- 
sity E, and gravitational constant G. Gravitational instabil- 
ity sets in when Q < 1. As magnetic fields are present in the 
discs, a magnetic Toomre parameter ( |Kim &; Ostri ker 2001 ) 

„ _ « (c s 2 +«j) 1/2 

Teg (5) 



can be defined as well, where va is the Alfvenic velocity tak- 
ing into account all components of the magnetic field. Since 
we know from Section [3.1| that the discs have almost Kep- 
lerian rotation velocities, we can easily replace the epicyclic 
frequency k in the above equations by the rotation frequency 
of the gas. In Fig. [II] we plot Q and Qmag for the run 2.6- 
Rot-MlOO at the same two times shown in Fig. [2] As can be 
seen, for both times Q and Qmag are above the threshold of 
1 indicating that the disc is indeed stable against gravita- 
tionally induced fragmentation. Only at r = 35 AU and t = 
15 kyr does Q drop below 1. However, since Qmag, remains 
clearly above 1 and since no disc fragmentation is recognis- 
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Figure 11. Toomre parameter Q (solid lines) and magnetic 
Toomre parameter Qmag (dashed lines) for the runs 2.6-Rot-M100 
and 2.6-Rot-M100-B at different times. 



able at that stage in the bottom right panel of Fig. [2] this 
indicates that the magnetic field has a stabilising effect on 
the di sc (see also|Hosking &; Whitworth|2004| Machida et al. 



2005 ] |Hennebelle Teyssier||2008| |Duffin fc Pudritz| [2009 
Seifr ied et al.|20lf ). We note that also for other times the re- 
sults shown in Fig.^^remain qualitatively unchanged. This 
supports our conclusion that the sink particles are formed 
outside the discs by core fragmentation as shown in Fig. [To] 
We emphasise that also for the runs 2.6-NoRot-M100, 
2.6-Rot-M100-C, 2.6-Rot-M100-p2, 2.6-NoRot-M300 and 
2.6-Rot-M1000, where more that one sink particle forms, 
fragmentation is not a consequence of a Toomre instability 
in the disc but due to the collapse of distinct gravitation- 
ally unstable regions. Hence, in these runs the occurrence of 
multiple protostellar systems is due to the capturing of by- 
passing protostars and not due to disc fragmentation. The 
only exemption is run 2.6-Rot-M100-B where in total 36 sink 
particles form. In this run after roughly 11 kyr disc frag- 
mentation occurs around the most massive sink particle at 
a radius of about 25 AU. In this case both the classical and 
the magnetic Toomre parameter clearly drop below 1 shortly 



before the fragment forms as shown in Fig. 11 (green lines) 



4 DISCUSSION 
4.1 Disc formation 

4-1.1 Magnetic field misalignment 



In Section [3.11 we have shown that the first disc-like struc- 
tures with noticeable rotational support build up around 
5 kyr and are already well developed around 7.5 kyr (see 
Figs. [l]-|6|. Since the mass-to-flux ratio /i, however, stays 
well below the critical value of 10 found in non-turbulent 
simulations over at least the first 10 kyr, no rotationally sup- 
ported discs should be expected to form. Thus, the build-up 
of Keplerian discs from the very beginning on has to be at- 
tributed rather to the disordered magnetic field structure 
and the lack of a coherent rotation structure in the sur- 
roundings of the discs. The required amount of angular mo- 
mentum is mainly carried by shear flows produced by the 
turbulent motions. Recently, Hull et al. (2012) reported a 



significant misalignment of outflows and magnetic fields in 
low-mass star forming regions. The misalignment was mea- 
sured on scales of ~ 1000 AU, thus on the same scale consid- 
ered in Fig. [7] in this work. As discussed in Section pT2] our 



simulations suggest that the observed misalignment could 
also be a consequence of an unresolved disordered magnetic 
field structure. Moreover, since misaligned magnetic fields 
are simply a particular form of a disordered field structure, 
we suggest that the formation of Keplerian discs in the mis- 



aligned case (Hennebelle & Ciardi 2009 Ciardi & Hennebelle 



2010||Joos et al.|2012) and in our work can be traced back 
to the same physical origin. 

We note, however, that very recently Li et al. (2013) 



tried to confirm the results found by the authors mentioned 
before with a different numerical scheme. One striking re- 
sult of their study is that misalignment does not allow for 
Keplerian disc formation for \i < 5, i.e. for typical observed 
values. Moreover, based on the work of Joos et al. (2012) 



and Hull et al. (2012), Krumholz et al. (2013) estimated the 



expected fraction of Keplerian discs to be as low as 10-50 
%. In the presence of turbulence, however, we showed that 
in essentially any case - even for fi < 5 - Keplerian discs 
can form. 



4-1.2 Magnetic flux loss 



Recently, Santos-Lima et al. (2012b) suggested that the loss 



of magnetic flux due to turbulent reconnection is responsible 
for the reduced magnetic braking efficiency. Magnetic flux 
loss alone, however, does not necessarily change the classi- 
cal picture of magnetic braking, i.e. a well-ordered magnetic 
field and coherent rotating as studied in numerous simula- 
tions neglecting turbulence. Hence, in this case the value of 
[i would be an appropriate quantity to assess the magnetic 
braking efficiency. However, in our simulations \i stays well 
below the critical value of ~ 10 during the formation phase 
of the Keplerian discs (see Fig. [7| indicating that flux loss 
alone cannot account for the observed build-up of Keplerian 
discs. Thus, other mechanisms like the aforementioned dis- 
ordered magnetic fields and turbulent motions in the close 
vicinity of the discs - which are natural outcomes when us- 
ing realistic initial conditions for star formation - have to 
account for this. 

In a second paper Santos-Lima et al. ( 2012a| argue that, 
when considering scales of the size of the disc (~ 100 AU), 
[i will increase, a fact we also observe in our simulations. 
As mentioned in Section [3. 2| the increase of \i could partly 
be due to accretion along the field lines. Moreover, we con- 
sider it as crucial to reduce the magnetic braking efficiency 
already on larger scales (by means of disordered magnetic 
fields but not weaker field strengths as indicated by the low 
value of \i) since otherwise there would be no angular mo- 
mentum left on the small scales (100 AU) to form the discs. 
This idea is also supported by the findings that for a well- 
ordered magnetic field the inwards transport of angular mo- 
mentum is balanced by the (negative) magnetic torque at 
scales up to 500 - 1000 AU ( [Seifried et al.|2011) . Neverthe- 
less, we suggest that the magnetic flux loss occurring subse- 
quently in the evolution of the discs helps to maintain them 
in a rotationally supported state, in particular in the later 
phase (t > 10 kyr). 

Finally, we note that we did not cover the parame- 
ter space in a uniform way. Since we have a (at least) 5- 
dimensional parameter space spanning the strength of tur- 
bulence, magnetic fields and rotation as well as the core mass 
and the density profile, even a limited coverage of the pa- 
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rameter space would result in a considerable amount of dif- 
ferent initial conditions. In combination with the very high 
resolution of 1.2 AU this leads to significant requirements 
of computational resources making a uniform coverage of 
the parameter space impossible. Nevertheless, even with the 
simulations presented here it becomes clear that turbulence- 
induced disc formation works for a wide range of molecular 
core masses and for different turbulence strengths. Further- 
more, we showed that the formation of discs does not depend 
on the presence of an overall core rotation although the ori- 
entation of the discs does. 



4.2 Fragmentation properties 

The degree of fragmentation in most of the runs with 100 
Mq is relatively low (~ 5 sinks) given a highly gravi- 
tationally unstable core with supersonic turbulence. How- 
ever, 



Girichidis et al. (2011) have shown that the fragmen- 
tation properties of a molecular cloud core with a density 
profile p oc r -1-5 as chosen in this work strongly depend 
on the realisation of the initial turbulence field. Indeed, in 
run 2.6-Rot-M100-B the core fragments much more heavily 
forming 36 sink particle in the first 15 kyr. The formation 
of Keplerian discs, however, is not affected by the different 
fragmentation behaviour. 

However, since we do not include radiative feedback 
from the protostars in our simulations, in particular in the 
heavily fragmenting runs 2.6-Rot-M100-B, 2.6-NoRot-M300 
and 2.6-Rot-M1000 we might significantly overestimate the 
number of fragments ( Krumholz et al.]|2007 2010). For the 
runs 2.6-NoRot-M100, 2.6-Rot-M100, 2.6-Rot-M100-C and 
2.6-Rot-M100-p2, however, we find a similar number of frag- 



ments as in comparable simulations of Krumholz et al. ( 2007 



2010) (5 - 10 sinks), despite the fact that we did not include 
radiative feedback. This is most likely a consequence of the 
presence of strong magnetic fields. This was already stated 
by |Myers et al. ( 2012 ), who find that magnetic fields seem to 



suppress fragmentation in the diffuse core regions whereas 
radiative feedback mainly acts in the dense disc region. 

However, also in our simulations without radiative feed- 
back we find that disc fragmentation is not the main frag- 
mentation mechanism in our simulations. In all runs the sink 
particles are formed via the gravitational collapse of distinct 



overdense condensations (Klessen et al.||2000 


Heitsch et al. 


2001 


|Padoan & Nordlund |2004; Hennebell 


e & Chabrier 


2008 


see also |Mac Low & Klessen||2004| for an overview), 



only in run 2.6-Rot-M100-B a few sinks are formed via disc 
fragmentation. Comparing Q and Q ma g shows that the mag- 
netic field can strongly contribute to the stability of the discs 
(Fig. 11 ). We note that whereas in related work the reduced 



degree of fragmentation was often simply due to the sup- 
pression of Keplerian discs formation |Hosking &; Whitworth 



2004 ||Machida et al |2005| |Hennebelle Teyssier|2008| |Duf- 



fin & Pudritz 



2009 



Seifried et al. 2011), here it is indeed 



a consequence of rising the Toomre parameter above the 
critical value for stability. Hence, our result is complemen- 



tary to_thework of|Bate| fl2009| ), |Offner et ah] ( |2009[ |2010 ) 

and Myers et al. (2012), who show that also radiative feed- 



back suppresses the fragmentation of protostellar discs. For 
the accretion discs around high-mass stars, however, [Peters 
et al. (2010 2011) have shown that ionization feedback and 



mentation, at least outside of radii of ~ 500 AU. Therefore 
it seems that the question of the importance of radiative 
feedback is not yet fully resolved. 

The capture and inwards migration of secondarily 
formed sinks in the high-mass cases presents a path to 
form tight binaries (separation ~ 1 AU). However, in run 
2.6-NoRot-M100 and 2.6-Rot-M100 the binary forms after 



roughly 5 kyr (see Fig. 10) with the lower-mass sink hav- 
ing a mass of ~ 0.1 Mq. Hence, by that time an extended 
first core might still exist (e.g. Masunag a et al.||l998 



Ma- 



sunaga &; Inutsuka 2000), which might have merged with the 
more massive sink to form a single protostar. We note that a 
more detailed study on the process of binary formation will 
be postponed to a subsequent study. 

4.3 Comparison to related work 

The initial conditions of run 2.6-NoRot-M300 were moti- 
vated by the recent work of Myers et al. (2012). In their 



radiative heating is insufficient to strongly suppress frag- 



work the authors study the collapse of a molecular cloud 
core with similar properties under the additional influence 
of radiative feedback. Using a spatial resolution of 1.25 AU 
these authors find a 40 AU-sized Keplerian disc. With our 
run 2.6-NoRot-M300, which has a spatial resolution similar 
to that of My ers et ah] ( |2012| ), we can confirm this result 
finding even two Keplerian discs. The first disc in run 2.6- 
NoRot-M300 shows heavy fragmentation and temporarily 
gets completely disrupted but quickly re-establishes its Ke- 
plerian structure after each disruption event. Furthermore, 
we note that Myers et al. ( |2012| ) do not find a Keplerian 
disc when they reduce the spatial resolution by a factor of 
8. We performed a resolution study for run 2.6-Rot-M100 
with a four times and eight times coarser resolution, i.e. 
dx = 4.7 and 9.4 AU, respectively. We find that with de- 
creasing resolution the discs become smaller and progres- 
sively sub-Keplerian. Although there is still an indication 
of a rotationally supported structure in the run with a 8 
times coarser resolution, the results of our resolution study 
strongly suggest that a spatial resolution of at least a few 
AU - depending on the actual numerical method used - is 
needed to properly simulate the formation of Keplerian discs 
in a strongly magnetised environment. 

Recently, the emergence of magnetic bubbles expanding 
radially outwards in the plane of the disc have been observed 
in non-turbulent simulations (e.g. Seifried et al.||2011 Zhao 
et al.||20lT] Krasnopolsky et aT]|2012 ). This feature, which 
occurs after a relatively short time (t < 10 kyr), is a conse- 
quence of the ongoing mass accretion and - due to the con- 
ditions of ideal MHD - a corresponding magnetic flux accu- 
mulation. Once the magnetic pressure in the centre is strong 
enough to overcome the gravitationally drag, a low-density, 
highly magnetised bubble forms thus reducing the flux in 
the centre of the disc. Clearly, this bubble-like feature does 
not show up in our runs including turbulence. Since the ac- 
cretion rates of the sink particles of 10 ~ 5 Mq yr _1 and a few 
10 -4 Mq yr _1 for the low-mass and high-mass runs, how- 
ever, are comparable to the non-turbulent runs, one could 
naively expect the same amount of magnetic flux to accumu- 
late in the disc centre eventually resulting in the appearance 
of a magnetic bubble. Since this is not the case, the mag- 
netic field has to diffuse outwards before it is accreted onto 
the sink. This could be achieved by processes like magnetic 
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reconnection and/or numerical diffusion in the already ex- 
isting, rotationally supported disc. This agrees with the fact 
that fi increases in time (see top panel of Fig. j\ and for 
smaller radii (see also Santos-Lima et al.|2012a|b ). A reason 
for an increased numerical diffusion of magnetic fields could 
be the lower infall velocities associated with Keplerian discs, 
which in turn give magnetic diffusion more time to operate. 
We note that the occurrence of numerical diffusion can also 
be motivated physically since in this regime non-ideal MHD 
effects like Ohmic dissipation or ambipolar diffusion are ex- 



pected to come into play (e.g.|Nakano et al.| 


2002| Mellon & 


Li|2009 Duffin & Pudritz|2009| |Dapp et al. 


|2011|). 



As pointed out in Section |3.4[ the masses of the discs 
range from about 0.05 Mq up to a few 0.1 Mq and thus are 
a factor of about 10 below the sink particle masses. A sim- 
ple density threshold criterion was used to define the disc. 
We note that recently |joos et ah] ( 2012 ) suggested a more 



sophisticated criterion including also the velocity structure 
to determine the gas belonging to the disc. However, since 
both criteria can only give an approximate mass for the disc, 
we decided to stay with our simple criterion^] For a more 
sophisticated comparison with observational data anyway 
synthetic observations produced with a radiation transfer 
code would we required. 

Nevertheless, the disc masses in our runs agree rea- 
sonably well with a number of observations ranging from 
low-mass (e. 



J0rgensen et al. 



2009; 



to high-mass protostellar objects (e.g. | Fuller et al. 



[Enoch et al. 


2011 


|Fuller et al. 


2001) 



However, in particular for the latter case discs with signifi- 
cantly higher masses up to a few Mq are observed, although 
these discs might correspond to structures without signifi- 
cant rotational support (e.g. Shepherd et al.||2001| |Chini 



et al | [20041 |Fernandez-L6pez et al.||2011| |Preibisch et al. 



2011| . Recently, Keplerian rotation in a Class protostellar 
disc has been observed up to a radius of 90 AU ( |Tobin et al.| 
2012) in accordance with our findings. In contrast, |Maury| 



et al. 



( 2010J) do not find Keplerian discs around Class 
objects. Hence, to unambiguously answer the question at 
which stage Keplerian discs form, high-resolution observa- 
tions e.g. with ALMA are required. Our work clearly shows 
that under realistic conditions, i.e. when including turbu- 
lent motions, early- type, rotationally supported protostellar 
discs can form in typical star forming regions. 



5 CONCLUSION 

This study was aimed at testing the turbulence-induced disc 



formation mechanism reported in Seifried et al. (2012) for a 
wider range of initial conditions. We performed several simu- 
lations of collapsing molecular cloud cores with a wide range 
of masses (2.6 - 1000 Mq) as well as different turbulence 
strengths ranging from sub- to supersonic. We showed that 
independently of the cores mass, the strength of turbulence, 
or the presence of global rotation in all cases Keplerian discs 
build up after a few kyr already. 

The discs typically have a diameter of 50 - 150 AU and 
masses of 0.05 up to a few 0.1 Mq. Interestingly, fragmenta- 
tion mostly occurs in distinct overdense regions of the core 



We checked that a more sophisticated criterion similar to that 



rather than in the discs due to the Toomre instability. We 
showed that magnetic fields significantly contribute to the 
stability of the discs. Tight binaries can still form by disc 
capturing and inwards migration of secondarily formed pro- 
tost ars. 

From the work on disc formation under the influence of 
magnetic fields done so far a common picture seems to arise: 
Due to the turbulent nature of star forming regions the clas- 
sical picture of a well-ordered magnetic field structure and 
coherent rotational motions, used by Mouschovias & Pale- 
|ologou| ( [1980| to assess the efficiency of magnetic braking, 
seems to break down. Under realistic conditions a disordered 



and/or misaligned magnetic field structure (see also [ Hen- 
nebelle fe Ciardi|2009||Ciardi fe Hennebelle|2010[|Joos et al. 



2012 



2013) as well as turbulent shear flows emerge, which 
cause the classical magnetic braking efficiency to drop and 
allow for the formation of Keplerian discs in essentially all 
cases considered here. The likelihood of disc formation in the 
turbulent case is therefore higher than that for purely mis- 
aligned fields in non-turbulent cores ( Krumholz et al.|| 2013 



of Joos et al. ( 2012| > does not change the results significantly. 



Li et al.|[2 013). During the formation of the discs the mass- 
to-flux ratio stays well below the critical value of 10, hence 
magnetic flux-loss alone cannot explain the formation of Ke- 
plerian discs. However, the subsequent increase of \i due to 
magnetic flux loss and mass accretion along magnetic field 
lines certainly helps to maintain a rotationally supported 
state during the further evolution of the discs. 
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